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Abstract 

Airway remodeling is a pathophysiologic process at the clinical, cellular, and molecular level relating to chronic obstructive 
airway diseases such as chronic obstructive pulmonary disease (COPD), asthma and mustard lung. These diseases are 
associated with the dysregulation of multiple molecular pathways in the airway cells. Little progress has so far been made in 
discovering the molecular causes of complex disease in a holistic systems manner. Therefore, pathway and network 
reconstruction is an essential part of a systems biology approach to solve this challenging problem. In this paper, multiple 
data sources were used to construct the molecular process of airway remodeling pathway in mustard lung as a model of 
airway disease. We first compiled a master list of genes that change with airway remodeling in the mustard lung disease and 
then reconstructed the pathway by generating and merging the protein-protein interaction and the gene regulatory 
networks. Experimental observations and literature mining were used to identify and validate the master list. The outcome 
of this paper can provide valuable information about closely related chronic obstructive airway diseases which are of great 
importance for biologists and their future research. Reconstructing the airway remodeling interactome provides a starting 
point and reference for the future experimental study of mustard lung, and further analysis and development of these maps 
will be critical to understanding airway diseases in patients. 

Citation: Najafi A, Masoudi-Nejad A, Ghanei M, Nourani M-R, Moeini A (2014} Pathway Reconstruction of Airway Remodeling in Chronic Lung Diseases: A Systems 
Biology Approach. PLoS ONE 9(6): e1 00094. doi:1 0.1 371/joumal.pone.OI 00094 

Editor: Daniel J. Tschumperlin, Mayo Clinic College of Medicine, United States of America 

Received November 18, 2013; Accepted May 22, 2014; Published June 30, 2014 

Copyright: © 2014 Najafi et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits 
unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. 

Funding: No current external funding sources for this study. 

Competing Interests: The authors have declared that no competing interests exist. 
* Email: amasoudin@ibb.ut.ac.ir 



Introduction 

Airway remodeling is a term used to describe the dynamic 
processes in obstructive airway diseases. It usually refers to 
epithelial layer injury followed by structural changes in the airways 
and lung architecture [1]. However, the cellular and molecular 
processes depend on the type and the state of disease and the 
patient. Consequences of airway remodeling could include a 
decrease in pulmonary function and reduced responsiveness to 
bronchodilator therapy. Airway remodeling is reported in complex 
diseases such as asthma, chronic obstructive pulmonary disease 
(COPD), and Mustard Lung as the main respiratory clinical sign. 
Also, progressive dyspnea and airflow limitations, mucostasis and 
mucosal inflammatory reaction are usually associated with airway 
remodeling [2,3]. Mustard lung has an irreversible pattern of 
airway obstruction like COPD [4] without any evidence of 
emphysema. It is resistant to anti-asthma therapy and an 
irreversible pattern of obstruction. Based on these similarities with 
asthma and COPD, mustard lung can also be a good model for 
evaluation of airway remodeling. 

There is a need for a holistic approach to decode the massive 
amount of data generated with modern biological approaches. 
Systems biology can integrate multilevel views of cell physiology 
data generated by low and high-throughput techniques into a 
comprehensive understanding of nonlinear molecular properties. 
Generation of high-throughput omics data, including genomics, 
proteomics and metabolomics enable us to simultaneously 



measure and analyze cellular components at any given condition. 
Currently, large databases of heterogeneous biological data are 
available including gene expression profiles (microarray, EST, and 
SAGE), interaction data, and catalogs of gene or protein functions. 
Also, many computational tools and algorithms have been 
developed to identify biological modules or pathways in the 
context of biological molecular networks [5] . Consequently, the 
systems biology strategy may be able to identify and construct 
novel pathways, and as such, is an emerging biological tool of great 
interest [6]. 

Although individual components of this molecular interaction 
data have been studied for decades, the accumulation of huge 
datasets to create molecular networks is a topical advance in the 
field of molecular medicine [7,8]. Moreover, recent progresses in 
molecular biology have highlighted the necessity of a systems 
biology approach. So, reconstruction and disruption of biological 
networks and pathways, including metabolic pathways, protein- 
protein interaction networks (PPI), signal transduction pathways, 
and gene regulatory networks (GRN), has been a valuable tool in 
the abstraction of biological concepts [8]. Most studies in this field 
have focused on the reconstruction, analysis and modeling of 
intracellular and extracellular networks [9]. This approach 
becomes more important when applied to polygenic diseases for 
complex etiologies [10,11], while disease or abnormal pathways 
such as airway remodeling are given less consideration. Analysis of 
disease pathways has the potential to elucidate the molecular 
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mechanisms underlying disease progression and response to 
treatment. Accordingly, novel genes, proteins and pathways are 
reported in complex diseases such as cancer [12], Alzheimer 
disease [11], atherosclerosis [13], and Parkinson's disease [14], 
and these can be understood by utilizing PPI network models 
combined with gene expression data. 

In this study, we attempt to describe the process of airway 
remodeling pathway in mustard lung [15]. Interestingly, more 
than 45,000 of 100,000 Iranian exposed patients are suffering 
from the late effects of sulfur mustard (SM; 2-bis-chloroethyl- 
sulfide) after almost 25 years post-exposure [4]. The chemical 
warfare agent sulfur mustard as a potent alkylating agent is highly 
reactive vesicant that can cause airway epithelial injury. Recent 
studies on Iranians of about 20-25 years in age after exposure to 
SM have shown the most common late complications in 
descending order of frequency are found in the lungs, eyes, and 
skin [4] . Damage to the epithelium layer is known as a key factor 
driving airway remodeling. Airway remodeling is the greatest 
cause of long-term disability among patients with combat- 
exposure to SM gas [16-18]. COPD and mustard lung are 
similar in clinical symptoms and signs, but differ at the molecular 
level and interactions between them. 

Accordingly, we have prepared a master list of mustard lung 
related genes. PPI and GRN networks were then generated and 
topological analysis was performed on merged network. Then, key 
signaling paths in networks were identified and mapped manually. 
The first reconstruction of an airway remodeling pathway in the 
mustard lung model is presented. A schematic workflow of our 
work is shown in Figure 1 . 

Materials and Methods 

Static network modeling plays a central role in systems biology. 
Therefore, the reconstruction process is highly related to our 
knowledge base in characterizing mustard lung. However, owing 
to the large amount of data, we need appropriate computational 
tools for gathering, constructing, and analyzing generated 
networks and pathways. A schematic workflow of our work is 
shown in Figure 1 and is described in detail in the main text. 

Data collection 

Three chronic sulfur mustard exposed and three unexposed 
individuals as the control group were selected for gene expression 
profiling. Airway biopsy was used for microarray gene expression 
experiment. Differential expression results were achieved during 
our former work by one way ANOVA test that is equivalent to t- 
test for two groups, (a P. value <0.01 was accepted as statistically 
significant after using Bonferroni correction for multiple tests) [19]. 
Statistical significance was validated using a two-tailed t-test 
assuming unequal variance, whereby significance was achieved for 
p<0.01. Therefore, 122 differentially expressed genes were 
identified by microarray experiment. The raw and process data 
were deposited in the http:/ /lbb.ut.ac.ir/Download/LBBsoft/ 
Mustard-Lung-Miacroarray. Furthermore, related published gene 
expression results (real-time PCR) were analyzed for completing 
the microarray gene list. This investigation attempted for the first 
time to identify gene profiles using whole human genome 
microarray chips, so as to recognize potential new target molecules 
and pathways involved in the pathophysiology of SM-induced 
airway remodeling. For comparison, there have been several 
microarray studies on animal models involving analysis of rodent 
pulmonary tissue [20] and mouse skin [2 1] after exposure to SM. 

Extensive literature surveys were then used to identify and 
confirm factors (gene, protein and metabolite) involved in airway 



remodeling. Literature mining was performed using keywords 
including sulphur mustard, mustard gas, mustard lung, and airway 
remodeling in Pubmed and google scholar databases without time 
limitation. Chronic mustard lung studies were selected for further 
detailed review (File S5). Finally, after removing redundant 
reports, 50 genes were extracted and merged to the 122 genes 
identified with microarray experiment (172 genes). iHOP [22] is a 
Web-based tool and Agilent literature search is a plugin in 
Cytoscape [23] that were used to automate analysis of abstracts in 
search for gene names. Finally, data were manually curated. In all 
cases, the HGNC official gene symbol was used to identify genes 
and proteins. Gene set annotation enrichment analysis also were 
performed using DAVID web tools (http://david.abcc.ncifcrf.gov) 
that provides a set of functional annotation tools for the genes 
categorized into Gene Ontology (GO) terms. 

Ensemble network and pathway reconstruction 

The network and pathways were reconstructed based on the 
master gene list and the molecular interactions documented in 
related papers and on-line interaction databases. The Pathway 
Resource List (http://pathguide.org) is a meta-database that 
provides an overview of more than 300 web-accessible biological 
pathway and network databases [24] . Human-specific Interactions 
of protein-protein interaction (PPI) data were abstracted from 
Biomolecular Interaction Network Database (BIND) [25], Data- 
base of Interacting Proteins (DIP) [26], Mammalian Protein- 
Protein Interactions Database (MIPS) [27], Human Protein 
reference Database (HPRD) [28], and Biological General Repos- 
itory for Interaction Datasets (BioGRID) [29]. Furthermore, Gene 
regulatory network (GRN) data were obtained from analysis of 
related microarray experiments and search in GRN database such 
as GeneMANIA [30] and Search Tool for the Retrieval of 
Interacting Genes/Proteins (STRING). To achieve the GRN and 
PPI networks, we entered the mustard lung gene list to these 
databases and got several GRN and PPI networks for our master 
gene list. Some of these databases don't accept a gene list; in these 
cases, o each gene was entered into the database manually, and 
resulting interactions added to the networks using Cytoscape 
software. Also, Cytoscape plugins were used for extracting, 
merging, visualizing, and analyzing unified interactive data. 
Molecular species in generated networks are represented as nodes 
and the interactions between these nodes as edges. 

Some of the achieved interactions were checked by pathway 
databases including Protein ANalysis THrough Evolutionary 
Relationships (PANTHER), Reactome, the National Cancer 
Institute Pathway Interaction Database (NCI-PID;), the Kyoto 
Encyclopedia of Genes and Genomes (KEGG), and Pathway 
Common, as the most widely used databases (table 1). 

Also, generated networks were compiled in Simple Interaction 
Format (SIF) amenable to Cytoscape for further topological 
analysis. After merging networks, the topological and statistical 
significance of the network was calculated using the Network 
Analyzer plugin in Cytoscape. The degree (connectivity) of nodes, 
the number of hubs (highly connected nodes), the shortest path 
lengths between any two nodes, the mean path length (the average 
of the shortest path lengths) and the network diameter (the 
maximum of the shortest path lengths) in comparison with random 
networks (Erdos-Renyi and Barabasi-Albert models; generated by 
random networks, a plugin in Cytoscape) were analyzed. 

Manual Pathway curation provides the most reliable means of 
extracting information from the literature, databases, and exper- 
iments. Thus, based on generated networks (PPI network and 
GRN), functional modules, signaling pathway databases, and 
previous knowledge, the airway remodeling pathway map was 
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Figure 1. A schematic picture demonstrating the pathway reconstruction workflow. The workflow shows that the master gene list is 
provided from two sources (microarray gene expression and literature mining). Protein-Protein interaction network (PPi) and Gene Regulatory 
Network (GRN) are compiled using Cytoscape software and its plugins. Also, signaling pathway reconstruction is accomplished in CellDesigner 
software based on biological knowledge accumulated in related literatures and public databases (e.g. Gene Ontology, Uniprot, and NCBI geneRIFs) 
and related signaling pathways deposited in the signaling pathway databases such as Biocarta, Reactome, and NCI-PID. The final mustard signaling 
pathway was drawn manually using CellDesigner. 
doi:1 0.1 371 /journal.pone.01 00094.g001 



manually constructed by CellDesigner version 4.2 (www. 
celldesigner.org) in Systems Biology Markup Language (SBML) 
format, which provide graphical environment, pathway visualiza- 
tion, and navigation. It can be used for building diagrammatic- 
based biochemical networks represented as a process flow, with the 
possibility of running functional simulations. 



Based on functional modules derived from mustard lung 
network, we used signaling pathway databases (Table 1) to map 
the functional modules to signaling pathways. Each path was 
identified in different signaling pathways and was manually 
curated and reconstructed (drawn) using CellDesigner. This stage 
requires multidisciplinary skills and a priori knowledge. 
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Table 1. Pathway and interactome databases used to identify interaction (edge) between genes in the master list. 



Database 


Version/date 


Graph type 


Address 


STRING 


V8.2 


Directed/Undirected 


http://string-db.org 


Reactome 


v43 


Directed 


www.reactome.org 


KEGG 


2012 


Directed 


www.genome.jp/kegg 


GeneMania 


3.1.1 


Undirected 


www.genemania.org 


NCI-PID 


2012 


Directed 


http://pid.nci.nih.gov 


Biocarta 


2012 


Directed/Undirected 


http://www.biocarta.com/ 


GeneGO 


V2.5 


Directed 


www.genego.com 


Pathway Commons 


v 3.1.17288 


Directed 


http://www.pathwaycommons.org 


PANTHER 


2012 


Directed 


www.pantherdb.org 



doi:1 0.1 371 /journal.pone.01 00094.t001 



Results 

Candidate gene list 

Following the workflow described in the Methods section for 
data retrieval, 172 genes (122 genes from microarray and 50 genes 
from literature) were selected according to databases, literature 
mining and related microarray data. The genes that are known to 
be involved in airway remodeling in sulfur mustard exposed 
patients are given in the File SI. These genes are annotated and 
described based on molecular process and functionality in the gene 
ontology (GO) and other annotation databases. 

Mustard lung networks 

We compiled a gene list (nodes) involved in airway remodeling 
based on PPI resources and literature mining. Briefly, PPI 
networks are commonly represented in an undirected graph 
format, with nodes corresponding to proteins and edges 
corresponding to physical protein-protein interactions. The PPI 
network obtained from the master list are given in Cytoscape .sif 
format in File S2. Based on current knowledge of physical 
interactions in PPI databases, we could identify interactions for 
just 96 of 172 nodes in the master list. Hence, a PPI network of 96 
nodes and 211 edges was generated. 

Gene expression data and co-expressed genes have become a 
useful resource in representing of the molecular state of cells. 
GRNs are commonly shown in a directed graph. GRNs consist of 
nodes, representing the expression profile of a particular gene, and 
edges representing significant associations between expression 
profiles (for example Pearson correlation). By using our gene list 
and co-expressed genes (extracted from GRN databases and 
related microarray data), a mustard lung GRN was generated. 
This network with 171 nodes and 1002 edges is shown in 
Cytoscape .sif format in File S3. Also the pdf files of networks are 
presented in File S6 (PPI) and File S7 (GRN). 

The mustard lung network was created by merging networks. 
The complex network is likely composed of several sub-networks 
or functional modules contributing to various diverse biological 
processes in mustard lung disease. Three major functional modules 
are illustrated in figure 2. 

Proposed signaling pathway 

Pathway extraction or static modeling of networks (PPI network 
and GRN) has become a very active area of research. Figure 3 
shows a schematic pathway that is generated by CellDesigner in 
SBML format (File S4). This pathway has been extracted from the 
generated mustard lung network. There are several signaling paths 



such as ERK/MAPK, EGFR activation, and EPAS1/ARNT 
paths in this disease pathway that have high crosstalk together. 

Network topology analysis 

Network parameters were performed by using NetworkAnalyzer 
plugin of Cytoscape for the mustard lung network. Topological 
analysis of a network identifies the global qualitative properties of 
the system. Network topology is used to provide the significance of 
a node in communicating with other nodes. Scale-free networks 
share two important functional characteristics. First, they are 
differently sensitive to damage. So if a small, peripheral node stops 
functioning, the network is very likely to continue working without 
problem. By contrast, if a hub is damaged, the functionality of the 
entire network is likely to be jeopardized. These topological 
characteristics are seen in biological networks [31]. In other words, 
if a hub node such as MMP9 is closed, most of the nodes and edges 
will be affected. Since biological networks satisfy power-law degree 
distribution, we have checked the mustard lung network for basic 
network properties. Degree distribution of a scale-free network 
having k connections to other nodes satisfies the following relation 
[32]: 

P(k)~k y 

where y is power-law parameter. For mustard lung network, we 
have applied curve fitting to aforementioned relation and 
calculated values of y and r 2 (coefficient of determination or R- 
squared). According to the r, the results show that mustard lung 
network is scale-free (Figure 4-right). The betweenness centrality 
Cb(n) of a node n is computed as follows: 



Cb(n)= ^^^^stW/ost) 



where s and t are nodes in the network different from n, o st denotes 
the number of shortest paths from s to /, and a st (rij is the number of 
shortest paths from s to t that n lies on (Figure 4-left). 

Topological parameters such as the average clustering coeffi- 
cient, topological coefficient, average degree and power law 
distribution of degrees and betweenness centrality have been 
assessed to predict the topology of the networks. The distribution 
of clustering coefficient is an important characteristic of biological 
scale-free networks. The clustering coefficient C„ of a node n is 
defined as C„ = 2eJ (kjk n -l)), where k n is the number of neighbors 
of n and e n is the number of connected pairs between all neighbors 
of n. Therefore, the clustering coefficient, characteristic path 
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Figure 2. Major disease-risk modules in mustard lung network. Colored nodes are selected based on fold change ranking in microarray gene 
list and previous reports. They have a central role in these modules as driver nodes. The A, B, and C sub-networks extract from the mustard lung 
network using Hubba, a plugin in Cytoscape software. 
doi:1 0.1 371 /journal.pone.01 00094.g002 



length, network centralization and network density of the Mustard 
lung network were compared with simulated randomized model 
networks (Erdos-Renyi and Barabasi-Albert models; generated by 
random networks, a plugin in Cytoscape) as shown in Table 2. 
The general structure of the network is far from random as shown 
by the comparison with the simulated network. In particular, the 
clustering coefficient of the Mustard lung network greatly differs 
from random network. These data indicate that the interactome 
possesses its own peculiar structure, and not a random one, which 
is a sign of a highly organized architecture. 

The clustering coefficient (table 2) is the average of the 
clustering coefficient for all nodes in each network. The 
characteristic path length gives the expected distance between 
two connected nodes in each network and the network density and 
network centralization show how densely the network is populated 
with edges and reflect the tendency of a network to contain hub 
nodes. 

Discussion 

Cellular signaling is a complex system governing the events in 
the organism that can be considered as communication both 
within and between cells. Airway remodeling occurs in many 



chronic lung disorders. The bronchial epithelium, sub-epithelial 
myofibroblasts and ASM (airway smooth muscle) cells are the 
major cell types involved in tissue repair processes in the airway 
tissue. Highly dynamic interactions and crosstalk between these 
resident cells via signaling proteins through autocrine or paracrine 
mechanisms dictate how tissue repair should progress. Signaling 
proteins from different pathways may interact directly (e.g. by 
phosphorylation) or influence each other indirectly (e.g. via 
regulation of gene expression). One component may also act in 
more than one pathway. Such cross-talk events can result in 
unexpected behavior. 

For the first time, we reconstructed an enriched network and a 
pathway that comprehensively catalogs major signaling pathways 
in the field of airway remodeling associated with sulfur mustard 
exposure. Disease pathway construction is essential to decipher 
complex regulatory patterns in deranged cellular processes. 
Therefore, identifying disease networks and pathways is important 
for understanding the pathology of disease and improving clinical 
diagnosis and treatment. Based on network topology analysis in 
comparison with random network, the generated networks have a 
larger clustering coefficient, higher network density, higher 
network centralization, and shorter characteristic path length 
than generated random networks. All of these features reveal a 
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Figure 3. The proposed pathway of airway remodeling in mustard lung. Pathway nodes are indicated by color coding; Green: single protein; 
White: protein complex; Yellow box: gene; Bisque: receptor protein; Purple: biological process. The pathway illustrates several paths such as ERK/ 
MAPK, ERBB1/ERBB2 complex activation via TFF3, and EPAS1/ARNT transcription factor activation. Different paths were extracted from pathway 
databases (table 1) and then were curated and reconstructed manually in the CellDesigner graphical interface. 
doi:1 0.1 371 /journal.pone.01 00094.g003 



high degree of connectivity found in biological scale free networks. 
This type of network with high levels of global connectivity and 
high-degree hubs is often found in biological systems that respond 
to changing external conditions [33]. However, while this system 
design is useful for normal tissue, in abnormal conditions such as 
mustard lung, the effects of network disruptions can spread 



rapidly. We have assembled a comprehensive map of putative 
regulators of airway remodeling in human airway epithelial layer 
of mustard lung. This pathway is composed of complexes of 
physically interacting genes and proteins, which are investigated 
with molecular biology approaches in genomics, proteomics, and 
metabolomics. 
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Figure 4. Analysis of mustard lung network with 1 72 nodes and 11 69 edges. These graphs were generated using NetworkAnalyzer plugin 
of Cytoscape. The scatter plot of betweenness centrality vs. number of neighbors (left) indicates that a limited number of nodes control the 
information flow between other nodes within the biological network. This means that a limited number of nodes with high interactions (hubs) 
control other nodes with lower interactions. The node degree distribution (right) shows that the network is scale-free considering the power-law 
degree distribution P(k)~k y (fitting result is y = 0. 790 and R-squared = 0.667). This means that the mustard lung network is a biological network which 
differs from random networks. 
doi:1 0.1 371 /journal.pone.01 00094.g004 



Previous studies have shown various signaling processes such as: 
angiogenesis, inflammation, oxidative stress, mucous secretion, 
apoptosis/ antiapoptosis, and ECM degradation in airway remod- 
eling [15,34,35]. Clinical and pathological signs of mustard lung 
also confirmed these processes. Disruption of a pathway may be a 
good thing. Thus, the detailed understanding of molecular 
pathways is an important step in diagnosis and treatment of 
disease. Eventually, it can help to determine biomarkers for 
complex diseases, which may lead to new approach for treatment 
and diagnosis [36,37]. 

Several signaling paths were identified in the generated mustard 
lung pathway, including EPAS1/ARNT path, ERK/MAPK path, 
MMPs path, and EGFR path. These highly connected paths can 
lead to various processes and phenotypes. Details and functions of 
these processes are discussed below. 

Extracellular Matrix and adhesion molecules regulation 

The extracellular matrix (ECM) is a complex composition of 
different combinations of collagens, proteoglycans, hyaluronic 
acid, laminin, fibronectin, and many other glycoproteins, includ- 
ing proteolytic enzymes involved in the degradation and 
remodeling of the extracellular matrix. One of the major roles of 
ECM is to serve as a structural framework for cell attachment and 
migration. Among proteases, matrix metalloproteinases (MMPs) 
are extracellular degrading enzymes that play a critical role in 
degrading components of the ECM and epithelial cell junctions, 
cancer metastasis, and liberation of tethered growth and chemo- 
tactic factors in response to wounding [38]. One of the lung 



disorders that has both inflammatory and remodeling features is 
mustard lung. Imbalance between levels of MMPs and TIMPs in a 
particular location can lead to accumulation of ECM protein at 
that site. The most well-known profibrotic cytokine TGF-(3 is 
known to crosstalk with the ECM pathway. Also, TGF-P signaling 
has a critical role in activation of angiogenesis [39] . 

Mucus hyper-secretion 

Mucus hyper-secretion is a major pathophysiologic feature in 
chronic inflammatory airway diseases. Trefoil Factor3 (TFF3) is 
one starting point in the airway remodeling pathway as shown in 
figure 3. Trefoil family has three homologous (TFF1, TFF2, and 
TFF3) [40,41]. Although the exact function of TFF3 in the 
respiratory tract remains largely unknown, experiments using 
recombinant human TFFs (hTFF) show that they may serve 
several roles in maintaining intestinal homeostasis such as by 
potentiating intestinal permeability, increasing the mucous integ- 
rity through mucin interactions, or directing epithelial restitution 
by promoting cellular migration. Also, Greeley et al. have shown 
that TFF1, TFF2, and TFF3 have roles in proliferation, as well as 
migration, cell death and differentiation phases of airway epithelial 
repair [42]. In our pathway, TFF3 acts as a ligand and binds to 
ErbBl/ErbB2 complex (ERBB1 also known as EGFR). Hetero- 
dimerization of ErbBl/ErbB2 stimulate ERK kinase cascade. 
ErbBl and ErbB2 recruit Son of Sevenless homolog (SOS) via 
adaptor protein GRB2 and SHC transforms protein (She), 
respectively. SOS is a guanine-nucleotide exchange factor for 
small GTPases, including H-Ras. Following this activation, the 



Table 2. Basic network parameters of the two generated networks, compared with simulated randomized model networks. 







Mustard lung network 


Simulated Barabasi Albert model, (scale free) 


Simulated Erdos-Renyi model 


Number of nodes 


172 


172 


172 


Clustering coefficient 


0.208 


0.127 


0.054 


Characteristic path length 


2.419 


3.206 


2.610 


Network density 


0.079 


0.023 


0.051 


Network centralization 


0.150 


0.225 


0.049 



doi:1 0.1 371 /joumal.pone.01 00094.t002 
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MAPK/ERK path is activated and we have shown MAP2K2 
over-expression that has a critical role in ERK pathway. Then, 
this cascade continues with the activation of downstream 
pathways. At the end of the pathway, Muc5B gene is expressed 
which leads to enhanced mucus secretion. Also, this cascade causes 
expression of TFF3 gene, leading to increased TFF3 protein in 
airway tissue. 

Angiogenesis 

Angiogenesis is a hallmark in the pathology of many diseases, 
including cancer, ischemia, atherosclerosis, and inflammatory 
diseases. Angiogenesis has additional roles in normal development 
and physiological processes in adults, including wound healing and 
tissue regeneration. Most previous studies have shown that 
hypoxia inducible factor 1 alpha (HIF1A) and endothelial PAS 
domain protein 1 (EPAS1) transcription factors promote the 
transcription of vascular endothelial growth factor (VEGF) when 
they interact with aryl hydrocarbon receptor nuclear translocator 
(ARNT). Moreover, in-vitro experiments showed that EPAS1 (or 
HIF2A) mRNA is expressed in the epithelial cells of pulmonary 
alveoli [43]. This complex binds to the promoters and then 
activates target genes involved in proliferation and differentiation, 
such as MMP9, VEGFA, CYP1B1, CCND1, and MUC5B. The 
two subunits EPAS 1 / ARNT are together forming a transcription 
factor that regulates expression of VEGF that enhances angio- 
genesis in airway tissue. Also EPAS 1 /ARNT regulates more than 
100 genes which regulate important mechanisms such as 
angiogenesis, apoptosis, and anaerobic metabolism. In mucus 
hyper-secretion and hypoxia, the EPAS 1 /ARNT complex is 
activated, and thus expression of downstream genes is increased as 
a result. EPAS 1 protein increases its own promoter activity as well 
as the transcription of the VEGF gene. 

It has also been noted that in prolonged hypoxia in lung 
epithelial cells, the up-regulation of HIFloc was transient, whereas 
increases in HIF2a were sustained [44] . A transient up-regulation 
of HIFlot could not be ruled out in these conditions, but the 
expression of HIFla was significantly lower than HIF2ot under all 
conditions. 

CXCL17 is a chemokine involved in tumor angiogenesis and is 
an anti-inflammatory factor that we have also observed to be over- 
expressed in mustard lung. CXCL17 induces the production of 
proangiogenic factors such as VEGFA. Also, CXCL17 is 
expressed in some aggressive types of gastrointestinal, lung and 
breast cancer cells [45]. 

TGF-betal is involved in extracellular matrix regulation, 
leading to airway remodeling in chronic airway diseases. Will- 
ems- Widyastuti et al. proved that TGF-betal acting through ECM 
regulation could also play a critical role in bronchial angiogenesis 
and vascular remodeling via VEGF pathway in asthma [46]. 

Inflammation and oxidative stress 

Oxidative stress plays a pivotal role in the pathogenesis of many 
chronic inflammatory lung diseases, especially in mustard lung. 
Oxidative stress induced cellular and molecular changes in the 
airway epithelium layer may significandy contribute to the 
pathogenesis of chronic inflammatory airway disorders and induce 
apoptotic or inflammatory responses in the lower respiratory tract 
[47]. 

Airway cells are constantly challenged by environmental (e.g. 
xenobiotics, UV, and drugs) and endogenous stressors (e.g. 
reactive oxygen species (ROS), hydroproxides, carbonyls, and 
quinones). If unchecked, these stresses lead to airway and lung 
inflammation. FM02/FM03 complex is one of the key enzymes 
that catalyze ROS to 02. With decreasing expression of this 



complex, ROS increases and therefore, leads to oxidative stress 
and inflammation. In addition, ROS is a second messenger 
required for regulating MUC5B expression and mucus secretion 
via EGFR signaling path [48]. 

Apoptosis and anti-apoptosis regulation 

Apoptosis (programmed cell death) and necrosis are the two 
major types of cell death. The process of apoptosis is triggered by a 
diverse range of cellular signals. Either the extrinsic or intrinsic 
death pathways can lead to apoptosis. Extracellular signals include 
toxins, hormones, growth factors, cytokines, nitric oxide, heat, 
radiation, nutrient deprivation, viral infection, hypoxia, the 
binding of nuclear receptors by glucocorticoids, or increased 
intracellular calcium concentration can damage DNA or cause 
cellular stress, triggering the release of intracellular apoptotic 
signals. There are some genes in File S3 that are connected to 
apoptosis, including CRIP2, PACS-2, and IFI35. Over-expression 
of CRIP2 induces apoptosis through activation of caspases 3 and 9 
[49]. CRIP2 acts as a transcription repressor of the nuclear factor- 
KB-mediated proangiogenic cytokine expression [50]. TRAIL/ 
DR5 triggers PACS-2 to traffic Bim and Bax to lysosomes to 
release cathepsin B and induce apoptosis [5 1] and IFI35 is induced 
by IFN-gamma [52]. Some of the identified factors (such as 
TSPAN8, CXCL17, CRIP2, STARD10, CSNK2A1) and path- 
ways (i.e. apoptosis, TGFfJ, VEGF and ERK signaling) are known 
to participate in cancer as well as their identification here linked 
with airway remodeling in mustard lung. Our data suggest that 
cancer and airway remodeling share several overlapping processes. 

This work has a few limitations. First, airway biopsy sampling is 
an invasive technique (especially for normal group), and most of 
the chronic mustard patients have about 25-year difficulties in 
their life style and strongly refuse to be included in the study. 
Combined with various exclusion criteria including age, smoking, 
and lung performance, the study population was quite small, thus 
gene profiling of each biopsy sample was performed in duplicate 
for further accuracy. Second, our analysis necessarily includes data 
from mouse models of airway remodeling which may or may not 
be representative of human disease. 

Conclusions 

Our findings demonstrate that the analysis and modeling of 
complex biological networks are beyond the capabilities of existing 
computational techniques that are needed for this type of scientific 
endeavor. Therefore, manual curation and pathway reconstruc- 
tion was followed after computational results. Current diagnosis of 
the airway remodeling diseases, in particular mustard lung and 
COPD, is mainly based on a combination of lung function 
evaluation and an observation of symptoms. There is no sole and 
exact clinical or laboratory test available. This highlights the need 
for biomarkers to diagnose disease or identify disease phenotypes. 
Therefore, the identification of key candidate genes, and their 
roles in regulating pathways, shows the need for bringing systems 
biology to the clinic as a powerful new approach. TFF3, ERBB2, 
EPAS1, and COL7A1 are the candidate factors identified here as 
potentially centrally involved in the pathophysiology of airway 
remodeling in mustard lung. The regulation of these proteins may 
be potentially useful in the treatment of airway remodeling. 
Finally, a comprehensive understanding of biological pathways 
can aid in the development of drugs to target specific cellular 
mechanisms while avoiding unwanted side effects. 
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Supporting Information 

File SI List of Genes that are known to be involved in 
airway remodeling in sulfur mustard exposed patients 
and their propertises. 

(XLS) 

File S2 The PPI network obtained from the master gene 
list are given in Cytoscape .sif format. 

(sir) 

File S3 The GRN network with 171 nodes and 1002 
edges is shown in Cytoscape .sif format. 

(sir) 

File S4 The proposed mustard lung pathway that is 
generated by CellDesigner in SBML format. 

(XML) 
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